home *** CD-ROM | disk | FTP | other *** search
/ Shareware Grab Bag / Shareware Grab Bag.iso / 007 / modula.arc / AFAI.MOD < prev    next >
Encoding:
Text File  |  1985-05-30  |  1.1 KB  |  50 lines
  1. IMPLEMENTATION MODULE afai;
  2.  
  3. CONST eps = 1.0E-8;
  4.  
  5. PROCEDURE sqrt(x:REAL):REAL;
  6. VAR a,c: REAL;   (* 0<x<2 *)
  7.  
  8. BEGIN
  9.   a := x; c := 1.0 - x;  (* abs(c) < 1 *)
  10.   REPEAT  (* a^2=b*(1-c)>=(a*(1+c/2))^2=b*(1-c)*(1+c/2)^2=b*(1-0.75*(c^2)-0.25*c^3)) *)
  11.     a := a*(1.0 + 0.5*c);
  12.           (* a^2 = b*(1-0.75*(c^2)-0.25*(c^3)) = b*(1-(c^2)*0.75+0.25*c)*)
  13.     c := c*c * (0.75 + 0.25*c);
  14.           (* a^2 = b*(1-c) *)
  15.   UNTIL ABS(c) < eps;
  16.   RETURN a
  17. END sqrt;
  18.  
  19. PROCEDURE log(x:REAL):REAL;
  20. VAR a,b,sum: REAL;  (* 1<=x<2 *)
  21.  
  22. BEGIN
  23.   a := x; b := 1.0;
  24.   sum := 0.0;
  25.   REPEAT
  26.     (* log(x) = s + b*log(a), b <= 1, 1 <= a < 2 *)
  27.     a := a*a; b := 0.5*b;
  28.     IF a >= 2.0 THEN
  29.       sum := sum + b;
  30.       a := 0.5*a;
  31.     END
  32.   UNTIL ABS(b) < eps;
  33.   RETURN sum
  34. END log;
  35.  
  36. PROCEDURE recip(x:REAL):REAL;
  37. VAR a,c: REAL;       (* 0<x<2 *)
  38.  
  39. BEGIN
  40.   a := 1.0;
  41.   c := 1.0-x;
  42.   REPEAT             (* a*x = 1-c, abs(c) < 1 *)
  43.     a := a*(1.0+c);  (* x*a = (1-c)*(1+c) = 1 - c^2 *)
  44.     c := c*c;        (* x*a = 1-c *)
  45.   UNTIL ABS(c) < eps;
  46.   RETURN a           (* recip = 1/x *)
  47. END recip;
  48.  
  49. END afai.
  50.